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ABSTRACT 

This paper presents a method for estimating the forces and moments on the six 
axes of an Unmanned Underwater Vehicle (UUV) operating in close proximity to a host 
vessel. Forces and moments encountered by the vehicle in the host vessel launchway, near 
the vessel and clear of the launchway, and far away from the host vessel, are modeled 
from experimental and theoretical sources. A quasi - steady approach to the vehicle 
launch dynamics allows coupling of experimental data, potential flow theory, and 
hydrodynamic maneuvering coefficients to develop an algorithm capable of modeling 
vehicle and host vessel interactions with any user defined trajectory. 

A vehicle motion simulation study was performed using the Revised Standard 
Equations of Motion (EOM) augmented to incorporate the interaction force model for 
proximity effects of the UUV to a host. Application of the trajectory simulation and force 
and moment time history are developed for a underwater vehicle launching from a 
submerged host illustrating the hydrodynamic interaction of the vehicle to the host. The 
method is developed to account for host forward speed, UUV propulsor effects, and open 
water hydrodynamics. The robustness of the method is demonstrated for several launch 
trajectories. 

Thesis Supervisor: Dr. Jerome Milgram 
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Chapter 1 INTRODUCTION 



1.1 Motivation 

Underwater vehicles find applications in an increasing number of diverse missions 
and roles. Recent uses of AUV’s and ROV’s both in commercial and milit ary 
applications continue to outline the capabilities, and limitations of employing unmanned 
underwater vehicles (UUV). Integral to the use of a UUV, whether autonomous 
underwater vehicle (AUV) or remotely operated vehicle (ROV), is the process used for 
launch and recovery. This process is largely driven by the type of support ship, land based 
launch and recovery facility, or submerged host platform. Launch and recovery from a 
submerged host represents an opportunity to add flexibility to the vehicle mission, but 
adds complexity to vehicle design. The ability of a submerged host vessel, hereafter 
generically referred to as the host, to deploy a ROV or AUV to a mission commencement 
location in latitude, longitude and depth offers unique advantages in mission duration, real 
time interrogation, and covert mission accomplishment over other launch platforms. The 
ability to properly control a vehicle as it transitions from within the submarine launchway 
to a point far away from the host and back until finally rehoused represents a complex 
dynamic process for which improved models are necessary to investigate the forces, 
moments and trajectory of the vehicle during the course of its mission and through launch 
and recovery. The design of effective vehicle control is the goal of the hydrodynamic 
modeling of the launch and recovery process. 
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1.2 Research Objectives 



Launch and recovery of a vehicle from a submarine is expected to represent the 
most demanding case for control effectiveness. The effect of the host body speed and 
launchway geometry on the hydrodynamic flow characteristics in and around the 
launchway have been documented by experiment [11], Additionally, interaction effects 
between the UUV and the host have been quantified by experiment for specific cases [11], 
One can foresee that the combined effects which in turn are coupled, lead to a complex 
hydrodynamic modeling problem. 

The first objective of this thesis is to use experimental observations from a scale 
UUV within and in the vicinity of a launchway to develop a model to describe the forces 
and moments acting on a UUV as a function of hydrodynamic and geometric criteria 
relative to the host. Due to the specific nature of vehicle and host geometry on the forces 
experienced, the model will be tailored to the host and UUV. However this specificity is 
necessary only in so far as to proceed to the next objective, that of determining a 
reasonable simulation to determine the time dependent forces and moments acting upon 
the vehicle. 

This second objective is to incorporate this force model into a form of the Revised 
Standard Submarine Equations of Motion [5] to determine the forces and moments acting 
upon the vehicle given a specified trajectory. This investigation relies upon mating 
experimental data in the complex region within the launchway, with inviscid flow code 
results external to the launchway and in proximity to the host "boundary", and accounting 
for viscous effects of cross flow drag associated with the vehicle transitioning through a 
region where the vehicle axis is not aligned with the mean free stream velocity. The 
resulting control force simulator is intended to be adaptable to alternate vehicle 
geometries. 
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Chapter 2 FORCE MODEL 



2.1 Introduction 

The forces acting upon a UUV operating in the open ocean differ from those 
experienced when operating in the vicinity of a boundary, such as the free surface, ocean 
bottom, or host vessel For the purposes of this thesis, boundary interaction effects are 
limited to a plane rigid wall approximation of a host vessel except in the internal region of 
the launchway of the host. Interaction forces during launch and recovery, and through the 
transition to the "infinite fluid", are of direct interest in determining the feasibility of the 
control system and/or the operating envelope for launch and recovery. Significant work 
has been undertaken in the past, including computational modeling of these forces to 
assess the effects of high speed ejection of torpedoes from submarines [2], In the region 
of interest for UUV launch and recovery, it will be necessary to operate within a different 
regime of slow forward speed of the host and slow relative speed of the UUV to the host. 
This combination of slow host and UUV relative speed indicates the plausibility of 
modeling the UUV launch and recovery interaction forces as quasi - steady, the UUV 
stationary with respect to the host and slowly changing with respect to its variation in 
geometric relation to the host. This approach allows use of experimental results to 
develop a model for the forces in this region. Steady flow potential based codes may then 
be used to model the interaction effects of the vehicle near the host to resolve areas where 
experimental data is incomplete. 

2.1.1 Background 

The determination of the hydrodynamic forces present on a vehicle transitioning 
from a submarine launchway is complicated by the number of parameters involved and the 
complexity of the fluid -body interactions. For a baseline case with no vehicle in the 
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launchway Sjoblom and Schwemin [15] showed that the flow in a typical submarine 
launchway is unsteady and turbulent, with turbulence intensities within the cavity ranging 
from 10 - 20% of free stream velocity. Experimental observations document average flow 
magnitudes 0 - 10% of the free stream at the aft end of the cavity and approximately 30% 
of free stream at the midpoint of the launchway shutter. They document a large 
recirculation region in mid - region of the cavity and a separated region approximately 
30% of the cavity length forward. As reported, the turbulence shows a flow periodicity 
with frequency content between 1 and 4 Hz. External to the launchway and turbulent 
boundary layer potential flow predictions compare well with measurements. 

In order to assess the force and moment characteristics of a vehicle in launch or 
recovery transition, both the region within the launchway and external to it require 
modeling. The model must be sophisticated enough to capture three dimensional 
geometry effects, body - to - body interactions, as well as viscous flow effects. Inviscid 
potential flow models capture the proximity interaction effects; however, they do not 
include viscous drag effects including cross flow drag or lift forces generated by control 
surfaces and vorticity shed into the wake from separation. Numerical potential and 
viscous flow techniques can not at present handle all these effects so a combination of 
experimental results, numerical techniques, and semi-empirical non linear motion 
derivatives are suggested for investigating the launch and recovery envelope. 

2.1.2 Experimental Data 

Launchway Zone 

Jensen and Milgram [11] describe an experiment conducted using a 1/10 scale 
model of an UUV transitioning through a submarine launchway. Figure 2- 1 illustrates the 
basic configuration and geometry associated with the vehicle and the host. Principle 
dimensions of the UUV vehicle are listed in table 2- 1. 
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UUV Princip 


le Dimensions 


Length (ft) 


10.5 


diameter (ft) 


1.604 



Table 2-1 Principle Vehicle Dimensions 





Figure 2- 1 

Data was obtained for the forces and moments at various points of advance along the 
launchway centerline. Reference to a "guide can" structure existing in the torpedo tube 
launchway is made in the model tests. The vehicle was advanced from a position farthest 
aft, corresponding to the aft end of the propulsor duct flush with the forward edge of the 
guide can, to a position 10 inches forward of this reference (100 inches full scale). This 
range of traverse starts with the vehicle fully within the launchway and extends to a 
position with the aft end of the propulsor duct approximately even with the forward most 
point of the launchway. Data in all six degrees of freedom were obtained at each position 
along the launchway for different onset flow conditions corresponding to the speed of the 
free stream past the host allowing the effect of forward speed on the host to be modeled. 
The UUV model was configured with a ducted propulsor which was operated during the 
tests at various RPM for each advance position and onset flow combination to simulate 
the influence of propulsor forces and race effects on the induced flow. 
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Near Boundary Zone 

Once the vehicle is outside of the launchway, the geometry of the vehicle and the 
host can be approximated as a body in proximity to a plane rigid boundary. The approach 
taken here, utilizes results from a Boundary Integral Element Method (BIEM) panel code 
to determine the force and moment effects due to proximity operations external to the 
launchway. This method uses inviscid potential flow theory and does not include viscous 
effects. A description of the method and implementation is contained in [14], The effect 
of the propulsor on the induced flow, and accompanying modification to the pressure 
integral is not captured by the BIEM in this application. An attempt to include these 
effects is made through use of experimental data obtained for the UUV at various onset 
flow speeds, propeller RPM, and clearance from the approximated plane boundary. 

2.2 Modeling Results 

2.2.1 Coordinate System 

The force and moment models are referred to a right-hand orthogonal coordinate 
system of moving axes, fixed in the UUV with the origin located at the center of mass 
(CG) of the vehicle. The positive directions of the axes are x - forward, y - starboard, and 
z - downward in accordance with Standard Nomenclature [16]. 

2.2.2 Launchway Model 

The experimental data lead to a model of the forces and moments experienced by 
the vehicle as a function of specific dependent variables. The parameters varied during the 
testing used to develop the model were: 

(1) onset flow speed 

(2) advance along the launchway axis 

(3) propulsor rpm 
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(4) thruster flow and direction 



The data runs with the thruster running were not used in this analysis in order to limit the 
force and moment effects to those induced by body and launchway geometry and 
propulsor effects. The forces and moments may be developed in accordance with the 
following dimensional form: 



F = f(d,U,p,x,y,z) 

M = f(d,U,L,p,x,y,z) 

The following non-dimensional terms are defined for use in developing a model 
independent of scale. 



Advance 



Length 

F 

CFY = 

| PVof'A, 



CMZ = 



M 






A,L 



PTF = 



(PROPULSIVEFORCE ) 

-~pUof‘Aj 



where: CF = Force Coefficient 

CXI = Xfoment Coefficient 
PTF = Propulsive Thrust Factor 
Uqj = Onset flow speed 
Af = Frontal Area 
L = Vehicle Length 



For each data set the forces and moments recorded were converted to standard body 
coordinates with forces and moments limited to the horizontal plane being modeled in this 
case. A propulsion model was required which would allow robust modeling of propeller 
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effects. The non-dimensional Propulsive Thrust Factor (PTF) was used, with PTF 
determined for each data set by measuring the mean recorded thrust at zero RPM for each 
onset flow condition and advance to length position and subtracting this value from the 
axial force. The zero PTF data corresponds to the axial drag force experienced in the 
launchway due to the body position alone and was modeled as well. The procedure 
selected to arrive at a model recognized that significant variance existed in the 
experimental data. Figures 2.2 through 2. 14 reflect the non dimensional data with 
associated best least squares linear regression fit. 



Drag Force CFX vs A/L [(xpos-H)/l] 




Figure 2.2 Drag Force 
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CFY(PTF) CFY(PTF) 



CFY(PTF) for A/L 0.68 




PTF 

Figure 2.3 Lateral Force for A/L = 0.68 

CFY(FTF) for A/L 0.84 




PTF 

Figure 2.4 Lateral Force A/L = 0.84 
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CFY(PTF) CFY(PTF) 



CFY(PTF) for A/L 1.0 




Figure 2.5 Lateral Force A/L =1.0 

CFYfPTF) for A/L 1.16 




Figure 2.6 Lateral Force A/L =1.16 
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CFY(PTF) CFY(PTF) 



CFY(FTF) for A/L 1.32 




Figure 2.7 Lateral Force A/L = 1.32 

CFY(PTF) for A/L 1.48 




Figure 2.8 Lateral Force A/L = 1.48 
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CMZ(PTF) CMZ(PTF) 



CMZ(PTF) for A/L 0.68 




PTF 



Figure 2.9 Lateral Moment A/L = 0.68 

CMZ<FrF) for A/L 0.84 




Figure 2. 10 Lateral Moment A/L = 0.84 
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CMZ(PTF) CMZ(PTF) 



CMZ(PTF) for A/L 1.0 




Figure 2. 1 1 Lateral Moment A/L =1.0 

CMZ(PTF) for A/L 1.16 




Figure 2. 12 Lateral Moment A/L =1.16 
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CMZ(PTF) CMZ(PTF) 



CMZ(PTF) for A/L 1.32 




Figure 2. 13 Lateral Moment A/L = 1.32 

CMZ(PTF) for A/L 1.48 




Figure 2. 14 Lateral Moment A/L = 1.48 
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A One Way Analysis of Variance (ANOVA) was conducted on each data set of force and 
moment as a function of A/L and PTF. Table 2-2 summarizes the results. 



A/L 


Data 


F Ratio 


Prob > F 


r 2 


0.68 


CFY 


4.0022 


0.0627 


0.200087 




CMZ 


539.6004 


<0.0001 


0.971202 


0.84 


CFY 


2.6603 


0.1339 


0.210132 




CMZ 


2.2577 


0.1639 


0.184186 , 


1.0 f 


CFY 


276.3695 


<0.0001 


0.802538 




CMZ 


44.7450 


<0.0001 


0.396869 


1.16 


CFY 


88.9794 


<0.0001 


0.847589 




CMZ 


167.6809 


<0.0001 


0.912892 


1.32 


CFY 


166.9032 


<0.0001 


0.830764 ! 




CMZ 


116.0631 


<0.0001 


0.773429 


1.47 


CFY 


31.1761 


<0.0002 


0.757141 




CMZ 


857.6477 


<0.0001 


0.988475 



f (2) outlying points omitted 
Table 2-2 One Way ANOVA Results 

In Table 2-2 , the F-ratio is a statistical measure indicating whether the model is a 

better fit than the response mean, as determined by 

( 2> 2 / Degrees of freedom ) 



F - ratio= 



model 



/ Degrees of freedom) 

The category represents the amount of total variation accounted for by the model with 
the remaining variation attributed to random error and is determined as 

( 2 > 2 ). 



? - 






model 

Ctotal 



To conclude from Table 2-2 that the force and moment are dependent upon PTF the F- 
ratio must exceed an observed significance probability (P- Value) of 0.05. To determine 
the effectiveness of the model, is used in conjunction with the F-ratio. The case for 

A/L of 1.0 utilized the parameter to omit two outlying points from the model. The 
results indicate that for the majority of sample conditions the effect of PTF at each A/L is 
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present, and that a linear model is justified. Higher order models were not justified In 
order to arrive at a smooth model for force and moment as a function of PTF and A/L the 
dCFY(PTF)/ cPtf A f L is used to fit a 3rd order cubic function and the intercept is adjusted to 

make a best fit. Figure 2-15 and 2-16 reflect the best fit for force and moment within the 
launchway. It can be seen that a positive outboard force exists for a minimum A/L of 0.64 
and remains positive throughout the vehicle transition along the launchway centerline with 
a maximum value at an A/L of 1 . 15 which corresponds to approximately one quarter of 
the stem shadowed by the launchway. The lateral force then diminishes to 50% of this 
force at the maximum A/L of 1.48. Lateral moment effects are quite different, showing an 
inboard moment at minimum A/L followed by a zero crossing at A/L = 0.9. From this 
position the moment shows a continual increase to a maximum with the vehicle fully 
emerged from the host, corresponding in magnitude but opposite sign its value at 
minimum A/L. 



25 



aCMZ(PTF)/5PTF 



dCFY(PTF)/ePTF VS A/L 




A/L 

Figure 2. 15 Fit of xfy (ptf)/ cptf^ j l 

SCMZ(PTF)/dPTF VS A/L 




A/L 

Figure 2. 16 Fit of xmz ( ptf )J, cPtf a t ^ 
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The resulting force and moment model in the launchway are expressed as 

CFY = /j ( scfy ( ptf )/ cPTF a , l )a / L + f x {cFY Mercept ) 

CMZ = / 3 ( 3CMZ ( PTF ) j (PTF a ! l )aI L+f x [ CMZ Mercept ) 
where fj represents a linear fit 

^ represents a third order polynomial fit 



CFY (PTF) for A/L 1.0 




Figure 2. 17 Sample of Lateral Force Model Fit A/L =1.0 
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CMZ(FTF) for AIL 1.0 



1.0 

0.5 

0 

- 0.5 



- 1.0 

0 1 2 3 4 5 

FTF 

Figure 2.18 Sample of Lateral Moment Model Fit A/L =1.0 
2.2.2 Host Proximity 

With the vehicle near the boundary but not in the launchway, insufficient data 
existed in [1 1] to predict the effects of relative pitch and clearance from the boundary. 

To develop this model, an inviscid potential flow panel code was resorted to, in order to 
determine the near host interaction effects. Difficulties in the numerical hydrodynamics 
were encountered in modeling due to the blunt forward end of the UUV geometry. Initial 
paneled models attempted to reflect the bow accurately, however the large change in 
sectional area vs. length resulted in a large negative pressure at the bow leading to 
excessive axial suction force instead of the zero drag condition expected. Modification of 
the model to blend the change in bow cross section over increased length provided a 
model which corresponded well with strip theory and experimental data to a s imil ar 
vehicle [10], The resulting panel geometry is illustrated in Figure 2. 19. 



v 

V 



V 



V 




V 



oqn: y-( .153 * A/L 3 +.15e* A/L*+ .514* A/l~50g)*PTF, RK)J«, mu 4<v;LU 

Y linear Irasf S^ram Fit 
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Figure 2. 19 UUV Panel Geometry 

The proximity model is obtained in two parts. First, investigating the effect of relative 
proximity to the boundary both with and without relative pitch including far field 
determination of the Munk moment using potential flow. And second, investigating the 
effect of PTF on force and moment from experimental data. Figure 2.20 shows the 
suction force and moment as a function of gap/diameter, where gap refers to the distance 
from the vehicle molded shell to the boundary. The effect of pitch and gap/diameter on 
the force and moment is illustrated for selected cases in Figures 2.21 through 2.23. The 
non linear least square fit is superimposed in Figures 2.21 and Figure 2.22 showing good 
agreement over the range of combinations of pitch and gap/diameter. Figure 2.23 
represents the potential flow Munk moment determined with the vehicle in an unbounded 
fluid domain. 
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CFY (F/0.5*p*v *A) 



CFY and CMZ At Zero Rich 
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gap/dia 

Figure 2.20 Boundary Interaction At Zero Pitch 



CFY vs gap/dia for Pitch From 0° to 20° 




gap/dia 



Figure 2.21 Lateral Force Near a Boundary 
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CMZ vs gap/diameter for Pitch of T 
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Figure 2.22 Lateral Moment Near a Boundary 



Mirnk Moment Effect Far From Wall 




Pitch Angle [deg] 



Figure 2.23 Open Water Munk Moment 
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The effect of PTF on the force and moment is captured by data obtained in [1 1]. Figures 
2.24 and 2.25 represent the lateral force and moment effects due to a change in propulsive 
thrust factor. Significant scatter of experimental data was observed and several outlier 
points have been omitted in both curve fits, using the method described previously, to 
arrive at a model which makes physical sense and tends to the proper limit. 
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PTF 



Figure 2.24 PTF Effect on Lateral Force 
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Figure 2.25 PTF Effect on Lateral Moment 
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Chapter 3 SIMULATION MODEL 



3.1 Coordinate Systems 



The forces acting on an underwater vehicle are normally referred to a non-inertial 
coordinate system fixed in the body, however, motion of the vehicle is most easily 
described in an inertial reference frame. For the purposes of this analysis, the inertial 
reference frame is translating at the constant velocity of the host, assumed to be along the 
positive X axis, and therefore represents a host fixed reference. The Z axis is directed 
positive downwards and the positive Y axis is to starboard in this frame of reference. To 
obtain the forces and moments in the non-inertial reference it is necessary to transform the 
vehicle motion to this frame. The vector transformation between the reference frames is 
uniquely determined by three angles: <)) rotation about the X axis (roll), 0 rotation about 
the Y axis (pitch), vp rotation about the Z axis (yaw), which describe the angular 
orientation between the two frames and are termed the Euler angles. The order in which 
the transformation about these angular relations is observed is significant, and if taken in 
different sequence can produce a different vector result. In accordance with standard 
convention the transformation is accomplished by ordering the rotations as first yaw ( vp), 
then pitch (0), and finally roll(<J)). The transformation matrix describing the vehicle motion 
in translation in the non-inertial frame, given the inertial frame motion, is expressed by: 



T x (<!>,0, y/) = 



cos(#)cos(y/) cos(60sin( y/) -sin(60 

-sm(^)cos(^) + sm(^)sin(#)cos( y/) cos(^)cos( ^) + sin(^)sin(#)sin( y/) sin(0)cos(0) 
sin(^)sin( y/) + cos(^)cos( y/)sin(0) -sin(^) cos( y/) + cos($) sin(#) sin( y/) cos(^)cos(0) 



This transformation accounts for the rigid body rotation from host fixed reference to body 
fixed reference. The transformation of angular rate quantities follows the same Euler angle 
transformation convention and order of rotations. 
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The angular rate transformation is then: 






1 

0 

0 



0 

COS(0) 

-sin(^) 



- sin( 0) 
sin(^)cos($) 
cos(0)cos(0) 



The total transformation matrix may be represented as: 



T{<j>,9, y/) = 



0 

0 T 2 tf,0,y) 



With 9t 0 representing a vector in the inertial frame consisting of scalar quantities of u 
(surge), v (sway), and w (heave) velocities , and p (roll rate), q (pitch rate) and r (yaw 
rate), or their respective time derivatives; the previously defined Euler angles <{>, 0, and vj / 
are used to determine the vector representing the same vector 5R 0 , but in the body 
reference frame changed only in orientation. The relation 

X = T(t,8,y)9 i 0 



is used to determine the body relative velocity and acceleration in the body reference 
frame translating with the velocity of the host in the host axial direction. Finally, the body 
velocities are corrected for this constant translation of the body fixed reference frame to 
arrive at the total body relative velocities independent of the host. 



3.2 Equations of Motion and Implementation 



The equations of motion of a rigid body may be expressed for arbitrary origin in 
matrix form as follows: 
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This expresses the Revised Standard Equations of Motion [5] with the acceleration terms 
isolated to the left hand side. The right hand side consists of the sum of the gyroscopic 
and centripetal acceleration effects as well as forces arising from the fluid hydrodynamic 
effects on the vehicle. These last include forces and moments accounted for by 
maneuvering coefficients including control forces, hydrostatics, crossflow and vortex 
effects. The hydrodynamic forces are calculated in accordance with the standard 
equations of motion. These are described for the hydrodynamic terms as: 



X = X qq q 2 + X rr r 2 + X rp rp + X vr vr + X wq wq + X n v 2 + X ^ w 2 + X w iiw + X q uq + X uu u 2 

Y = Y P\P\ P \ P I + Y P1 Pq + YrUr + Y P UP + YvUV + y rH V M + Y ™ V2 + Y »P Wp + y vH V W + Y r\rf\ r \ + Y - UU + Y v U V 






■Jv 2 + 1 



Z = Z q uq + Z„uw + Z^u\w\ + Z vp vp + Z.uu + Z q ^q\q\ + Z vv | iv |H'|h-| + Z rf yp + Z Mv w 
K = K p up + K r ur + K qr qr + K^wp + K wr wr + K v uv + K p \ p \p\p\ + K.u 2 
M = M q uq + M q ^q\q\ + w| w\ + M vp vp + M rp r + M,u 2 + M w uw + M^, wyfv 2 + 
N = N pq + N up + N wp+N.u 2 + N r ur + N v uv + N<,r 



Vv 2 +W 2 



+ N , v ,u|v| + N w v 2 + jVi iv|v| + AT, ,r |r| 



The body terms, including coriohs, centripetal, and hydrostatic effects are calculated 
separately. Neglecting cross products of inertia these terms may be expressed in the form: 

X body = ~m(wq -vr - x G (q 2 + r 2 ) + y c pq + z G (pr))-(W - B)s\n(0) 

Y body=- m {ur-wp +x a pq + z a qr) + {W -B) cos(0)sm{<j>) 

Z body = ~m[vp - uq + x a pr - z G (p 2 +q 2 ) + (W - B)cos(0)cos(<fi) 

K body = (!yy ~ 1 zz )<l r ~ m (yc (vp - wq ) + z G (ur - wp )) + (y G fV -y g B) cos(0)cos( <t > ) - (z G W - z g B)cos(0)sin(0) 
M body =( I zz- I xx)P r ‘- m { x c( v P~ u< ])- z c( w ‘J ~ vr ))-( x G W ~ x b B)cos( 0) cos( <f>) - (z G (V - z g B) sin( 0) 
N body=( I xx - 1 yy)P0 ~ m ( x c( ur -™p)~ y G ( wt l -w )) + ( x G fv ~ x g B)cos(0) sin(^) + (y G fV -y B B)sin(0) 
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The hydrodynamic maneuvering coefficients may be obtained by experiment, or predicted 
by numerous methods [1], [9], [12], These relations do not include the propulsive thrust 
model forces or control forces in order to prevent redundancy of these terms with the 
interaction force model. 

The right hand side forces on the vehicle may be expressed as a partition of the 
total right hand side in the following fashion: 

Fioui = F H + Fi + Fc 

where F tolal may be equated to the acceleration and mass product terms , F H represents the 

forces and moments acting on the body due to open water hydrodynamic, centripetal, 
gyroscopic, and hydrostatic effects. The term Fj accounts for the interaction effects due to 
the operation in a restricted fluid e.g. in the launchway and host vicinity. F c accounts for 
the control force and moment required to effect the selected trajectory. F c may be solved 
for directly: 

Fc = F t0U | - (F„ + Fj ) 

This relation is adapted for implementation in a six degree of freedom motion simulation 
program for which a general description of the method employed follows. 

The desired trajectory of the vehicle must be known a - priori and specified as 
input to the algorithm. It is important in the launch and recovery of the vehicle to know 
the relative position , velocity, and acceleration of the vehicle with respect to the moving 
host, therefore the trajectory is specified in a host fixed inertial frame. The vehicle 
accelerations and velocities are determined by a coordinate transformation from the 
moving inertial frame to the body fixed frame with the effect of forward speed of the host 
accounted for in the transformation. The total force at each time step is determined from 
the vehicle accelerations by solving the product of the mass matrix and body relative 
acceleration vector directly. A word of explanation is warranted with respect to the use of 
a motion simulation program. For the aspect of this thesis, with the specification of a 
trajectory kn own explicitly over time along with acceleration, and velocity, there is no 
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need to invert the mass matrix and integrate over each time step. The forces and moments 
may be computed directly at each time step allowing the algorithm to run very fast without 
growing errors due to a integration scheme. The hydrodynamic forces and moments 
including the body terms are calculated at each time step without consideration to the host 
presence. A separate subroutine implements the interaction force model previously 
developed, where forces and moments when the vehicle is located within the launchway 
are determined completely from the launchway model and include hydrodynamic and 
interaction effects. Once the vehicle has transitioned to a position outside of the 
launchway the interaction effects alone are determined. This condition reflects only a 
portion of the forces acting upon the vehicle, therefore to arrive at the total right hand side 
force vector when the vehicle is outside of the launchway, the hydrodynamic derivative 
based forces are added to the interaction effects. 

In order to achieve a smooth transition between the experimental data representing 
the total force on the vehicle to the force calculated by the hydrodynamic derivative based 
force vector, it is necessary to employ a smoothing procedure to blend the differences 
between the two. In order to accomplish this, a zone structure is defined with the 
following description: 

Zone 1 Inside the host launchway F(x,y,z) = Hydrodynamic e + Interaction e 
Zone 2 Within the transition region F(x,y,z) = Hydrodynamic} + Interaction} 

Zone 3 Outside of the transition region F(x,y,z) = Hydrodynamics + Interaction e 
Where F(x,y,z) is the force calculated by the force model and subscripts 
e? j and s refer to experiment, interpolated, and coefficient based 
semi-empirical data respectively. 

Zone 2 is where the smoothing routine operates. This zone is defined by a hemisphere 
with the center located at the launchway exit defined by x c ,y c ,z c with the diameter of this 
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4 ) 



A geometry limitation file containing the position information with regard 
to launchway geometry, host geometry and host speed. 

5) A trajectory file which includes vehicle acceleration and velocity scaling 
parameters and trajectory modification commands versus time. 

A linear interpolation scheme is used to determine the body position, velocity and 
acceleration changes over time from the commanded trajectory inputs. A partial list of the 
subroutines implementing the trajectory, force model and coordinate transformation is 
included as Appendix A. 
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Chapter 4 SIMULATION RESULTS 



4.1 Trajectory Runs 

Four different trajectories were selected for evaluation of the simulation. A copy 
of the input files is included in Appendix B. The trajectories were selected to test the 
blending between the force and moment models between Zone 1 and Zone 3. The 
trajectories selected were not intended to reflect in all cases the trajectory that would be 
optimal or would represent the real vehicle desired course but would illustrate the model. 

4.2 Straight Launch 

The first trajectory simulated was a straight ahead motion of the vehicle. Using an 

exponential time function acceleration described by: 

x(t)=l( l-e -0 - 04 ' 2 ) 

(— 0 . 04 / 2 ) 

= (4*0.04 & J ) 

The vehicle is directed out of the launchway and away from the host at a constant heading 
angle. The trajectory in host reference is illustrated in Figure 4.1. Figure 4.2 through 4.7 
show the resulting accelerations, velocities, forces and moments for this maneuver in the 
vehicle reference frame. 
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Figure 4. 1 Vehicle Trajectory 
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Figure 4.2 Velocity Profile 
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Figure 4.3 Acceleration Profile 
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Figure 4.4 Hydrodynamic Forces in Open Water 
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Figure 4.5 Hydrodynamic Forces Including Vehicle Interactions 
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Figure 4.6 Total Force 
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Figure 4.7 Control Forces 

Although unrealistic as a desired maneuver due to the large sustained crossflow drag 
resulting from the constrained yaw angle, this trajectory is illustrative of the interaction 
effects and transition from a region where interaction effects play a role to the region in 
the far field where the interaction effects are small. As such it is useful as a baseline 
condition. The location of the vehicle relative to the zones previously defined is indicated 
in the figures, where pertinent, by a zone indicator line which takes a value of 1, 2, or 3 
respective to the appropriate zone. 

The interaction effects of the launchway and host vehicle are apparent in the figures. As 
the vehicle nose passes the launchway exit (occurring very early in the trajectory at ~ 4 
sec) the lateral force can be seen to decrease as the nose proceeds into the transition zone. 
The effect on the lateral moment indicates an almost impulsive effect turning the moment 
from a negative value when the vehicle is fully within the launchway to a large positive 
moment as the vehicle fully exits. The magnitude of this moment depends upon the 
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interaction due to the boundary as well as the effect of PTF on the accelerated flow 
between the vehicle and the interacting boundary. As the boundary recedes, the moment 
decreases to the open water Munk moment due to the imposed angle of attack of the 
vehicle in this trajectory simulation. 

4.3 Launch with Yaw Maneuver 

The second trajectory for simulation yaws the vehicle to a zero relative heading 
angle after the vehicle has cleared the launchway. This represents a more realistic 
maneuver intended to investigate the effect of maneuvering to remove crossflow. The 
vehicle is yawed to a zero heading angle at time {t = 25 sec} and proceeds along a path 
parallel to the host for the rest of the simulation. This trajectory is illustrated in Figure 
4.8. Figures 4.9 through 4. 14 illustrate the velocities, forces and moments. 
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Figure 4.8 Vehicle Trajectory 
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Figure 4.9 Velocity Profile 
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Figure 4. 10 Acceleration Profile 



Force [lb] Force [lb] 



Traj_2 Open Water Hydrodynamic Force 



80 r 



r 



80 



60 



40 



N 

Y 

X 

Zone indicator 



• 60 



-40 



20 



I 

1 

i 

t 

\ \ 

S \ 

V 



20 



-20 



J -20 



25 



50 

Time [sec] 



75 



100 



Figure 4. 1 1 Hydrodynamic Forces in Open Water 
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Figure 4. 12 Hydrodynamic Forces Including Vehicle Interactions 
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Figure 4. 13 Total Force 
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Figure 4. 14 Control Force 
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The interaction effects near the host are apparent in the lateral force and moment as they 
were in the straight launch trajectory. In addition, acceleration effects during the 
maneuver can be seen to affect the required control forces, specifically in this instance the 
contribution of N, to the lateral moment. 



4.4 Launch with Yaw and Deceleration Maneuver 

The third trajectory selected for simulation illustrates a maneuver intended to 
reflect a realistic vehicle maneuver. In this simulation the vehicle transitions out of the 
launchway and subsequently decelerates to the desired speed once clear of the host. 
Figure 4. 15 shows the trajectory history. Figures 4. 16 through 4.21 illustrate the 
velocities, forces and moments. 
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Figure 4. 15 Vehicle Trajectory 



50 



[ft/sec 2 ] [ft/sec] 




0 20 40 60 80 100 



Time [sec] 



Figure 4. 16 Velocity Profile 
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Figure 4. 17 Acceleration Profile 
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Figure 4. 18 Hydrodynamic Forces in Open Water 
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Figure 4. 19 Hydrodynamic Forces Including Vehicle Interactions 
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Figure 4.20 Total Force 
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Figure 4.21 Control Force 
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Again as in the previous trajectory the significance of the acceleration effects are evident 
in the axial and lateral body forces during deceleration. The lateral moment peak occurs 
just after the vehicle has cleared the launchway and experiences the maximum interaction 
effect of the planar boundary in proximity to the vehicle in this configuration. 



4.5 Straight Launch At Zero Host Speed 



The fourth trajectory selected for simulation illustrates a launch from a stationary 
host. In this simulation the vehicle transitions out of the launchway maintaining a straight 
course and speed without onset flow. Figure 4.22 shows the trajectory history. Figures 
4.23 through 4.28 illustrate the velocities, forces and moments. 
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Figure 4.22 Vehicle Trajectory 
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Figure 4.23 Velocity Profile 
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Figure 4.24 Acceleration Profile 



55 



Force [lb] Force [lb] 



Traj_4 Open Water Hydrodynamic Force 
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Figure 4.25 Hydrodynamic Forces in Open Water 
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Figure 4.26 Hydrodynamic Forces Including Vehicle Interactions 
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Figure 4.28 Control Force 
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It is of interest to note the differences between this launch and the previous launches, with 
forward host speed. In this launch the interaction effects are clearly seen; the lateral force 
changing sign from an outward force when the vehicle is mostly within the launchway to a 
inward suction force once clear. The lateral moment also changes sign from a bow in 
moment when the vehicle is mostly within the launchway to an outward moment due to 
boundary interaction and suction region acting on the aft portion of the vehicle. The axial 
force reflects the acceleration to constant velocity with corresponding constant axial 
hydrodynamic drag, while the lateral moment approaches asymptotically to a value 
dependent upon PTF. The effect of forward host speed has a marked effect on the force 
and moment on the vehicle, from commencement of launch with the vehicle fully within 
the launch way, to its transition to open water. With a 4 feet per second onset flow 
velocity due to the host forward speed, the lateral force does not achieve sign reversal, 
and an outward force dominates throughout the maneuver sequence. The effect on 
moment is pronounced as well, reflecting a twenty fold increase in peak moment over the 
zero forward speed case. 

The small perturbation of the lateral moment observed as the vehicle becomes fully 
emerged from the launchway is qualitatively similar to the moment profile results 
described in [2] for a high speed torpedo ejection from a submarine launchway. This 
comparison is not intended to validate the results here but is offered for completeness. 

Due to the disimilar nature of the ejection dynamics compared to vehicle swimout, a 
quantitative comparison was not conducted. 
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Chapter 5 CONCLUSION 



5.1 Summary 

From the comparison of trajectory runs it is clear that the force and moment model 
implemented into a six degree of freedom motion algorithm can solve for control forces 
versus time for any desired trajectory. The mathematical model successfully makes non 
dimensional and integrates experimental data with potential flow panel code results, and 
with a parametric based non linear mathematical model, to span the trajectory regime from 
launchway to open water. Host vessel speed effects are incorporated. The smoothing 
technique accounts for the inevitable mismatch between the experimentally derived 
interaction model to the parametric, coefficient based, hydrodynamic forces in open water 
and allows for geometrically similar models to be tested using the existing experimental 
results. 



Based upon trajectories with changing accelerations, acceleration effects add 
significant information to the force simulation. Vehicle motion gain evaluation is enabled 
through the trajectory input file. 

The magnitudes of the required control forces are broken down in such a way so 
as to enable the contributing interaction forces, open water hydrodynamic forces, or body 
mass forces, to be compared and their dependence upon trajectory evaluated. 

5.2 Recommendations for Future Work 

The control force output from the model represents, in a sense, the magnitude of 
control actuated force to maintain the selected trajectory. A shortcoming in the model 
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exists in that it ignores interaction effects between the control surface or thruster, the 
vehicle, and proximity boundary. The fidelity of the simulation would be improved by 
incorporating control effects on the uncontrolled vehicle interaction force model. 

For a vehicle with a different propulsor, a means of effectively adapting the thrust 
model and resulting thrust, moment and torque effects would be useful. 

The basis for a maneuvering workbench exists which can be used to assist in 
control appendage sizing, control system design, and free swimming motion simulation in 
one simulation environment. 

The accuracy of the simulation and its representation of the real world can only be 
determined by comparing its output to full scale tests. These would indicate those areas 
captured by the model as well as areas requiring improvement. Additional scale testing 
should be conducted to provide a more robust mathematical model, and should 
incorporate thruster, propeller, and launchway off - axis effects. 
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Appendix A. Source Code Listing for Trajectory Simulation and Force Model 

This appendix includes source code for the implementation of the trajectory, and force 
model. Also included is a portion of the Equation of Motion (EOM) program used to 
calculate the hydrodynamic open water forces. 



63 



Q ********************************************************************** 

SUBROUTINE TRJ 

£ ***************************** ********************^y'(Jg***************** 



IMPLICIT NONE 

REAL*8 CTRN(6,6),RHS(6),EDOT(6),E(6),PDOT(6) 

REAL* 8 P(6),DT,TA,PDOTO,PDDOTO 

REAL*8 DS4,DC4,DS5,DC5,DS6,DC6,DT5,DS5DC6,DS5DS6 

REAL*8 PDDOT(6),RHSTOT(6),RHSC(6) 

REAL*8 A(6,6),WORK(6,6),ACOF(6),BCOF(6) 

REAL* 8 ACP(3),BCP(3),RHSFM(6),RTD,DTR 

REAL* 8 XC, Y C,ZC,XMIN,XMAX, YMIN,YMAX,ZMIN,ZMAX,D2,HRAD 

REAL *8 VRAD,UHOST,DIA,RHOC,ATLMAX,VLEN 

REAL* 8 FX,FY,FZ,MX,MY,MZ, AEDOT 1 ,RHS 1,PTF 

REAL* 8 ATL,VEX(2),TEMP1,TEMP2 

REAL* 8 P1A,P1T,P2A,P2T,P3A,P3T 

REAL* 8 P4A,P4T,P5A,P5T,P6A,P6T 

REAL* 8 PlDOT,P2DOT,P3DOT 

REAL* 8 P4DOT,P5DOT,P6DOT 

REAL* 8 PlDDOT,P2DDOT,P3DDOT 

REAL* 8 P4DDOT,P5DDOT,P6DDOT 

INTEGER* 2 I,J,JJ,K,JJJ,II,KODE,JR,JJR,JP,JJP,JYA,JJYA 

INTEGER*2 JX,JJX,JY,JJY,JZ,JJZ 

CHARACTER* 15 LABEL 

REAL* 8 DD(1800) 

COMMON /COMDD/ DD 

COMMON/FGEOM/ XC,YC,ZC,XMIN,XMAX,YMIN,YMAX,ZMIN,ZMAX, 
1 D2,HRAD,VRAD,UHOST,DIA,RHOC,ATLMAX,VLEN 
1 ACP,ACOF,BCP,BCOF 

SAVE J,JJ,JR,JJR,JP,JJP,JYA,JJYA,JX,JJX,JY,JJY,JZ,JJZ 
DATA J,JJ,JR,JJR,JP,JJP,JYA,JJYA /9, 9,2, 0,2, 0,2,0/ 

DATA JX,JJX,JY,JJY,JZ,JJZ /2, 0,2, 0,2,0/ 

EQUIVALENCE (DD( 1), DT) 

EQUIVALENCE (DD( 11), JJJ),(DD( 4), RTD) 

EQUIVALENCE (DD( 10), TA),(DD( 5), DTR) 

EQUIVALENCE (DD( 31),EDOT(l)) 

EQUIVALENCE (DD( 55), E(l)) 

EQUIVALENCE (DD( 61), P(l)) 

EQUIVALENCE (DD( 71), RHS(l)) 

EQUIVALENCE (DD( 77),PDOT(l)) 

EQUIVALENCE (DD(107),A(1,1)) 

EQUIVALENCE (DD(83),PDDOT(l)) 
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DIMENSION P 1 A( 1 0),P LT( 10),P2A( 10),P2T( 10),P3A( 10),P3T( 10) 
DIMENSION P4A( 1 0),P4T( 1 0),P5 A( 1 0),P5T( 1 0),P6 A( 1 0),P6T( 1 0) 

IF (TA EQ. 0.) THEN 

C READ IN THE GEOMETRY CONSTRAINTS 

WRITE(6,*) ' OPENING AND READING GEOMETRY LIMITATION FILE' 
OPEN( 1 1, FILE- GEOM.LEVT, STATUS-OLD') 

READ(11,*) XC,YC,ZC 
READ( 1 1,*) XMIN,YMIN,ZMIN 
READ( 1 1,*) XMAX,YMAX,ZMAX 
READ(11,*) D2,HRAD,VRAD,VLEN 
READ(11,*) UHOST,DIA.RIIOC,ATLMAX 
CLOSE (11) 

C READ IN FORCE MODEL COEFFICIENTS 

WRITE(6,*) ' OPENING AND READING LEAST SQUARES COEFFICIENT 
FILE' 

OPEN( 12,FILE='Z2COEF.EXF, STATUS=’OLD') 

READ (12,*) (ACP(H), U = 1,3) 

READ (12,*) (BCP(H), II = 1,3) 

READ (12,*) (ACOF(II), II = 1,6) 

READ (12,*) (BCOF(H), H = 1,6) 

CLOSE (12) 

C READ IN TRAJECTORY MODIFICATIONS 

WRITE(6,*) ' OPENING AND READING TRAJECTORY COMMAND FILE' 
OPEN( 14,FILE='TRAJ.CM, STATUS^'OLD') 

READ (14,*) (VEX(U), U = 1,2) 

READ (14,*) 

READ (14,*) 

DO 201= 1,10 

READ(14,*)P1A(I),P1T(I),P2A(I),P2T(I),P3A(I),P3T(I) 

20 CONTINUE 
READ(14,*) 

READ(14,*) 

DO 301= 1,10 

READ( 14,*) P4A(I),P4T(I),P5A(I),P5T(I),P6A(I),P6T(I) 

P4A(I) = P4A(I) * DTR 
P5A(I) = P5A(I) * DTR 
P6A(I) = P6A(I) * DTR 
30 CONTINUE 
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CLOSE (14) 

C COMPLETE MODEL ENTRY 
ELSE 
END IF 

C USE STANDARD INITIAL VELOCITY PROFILE 
C MODIFY TRAJECTORY BY SCALING OF U,V, AND W RELATIVE 
C TO THE HOST FOLLOWED BY ROLL,PITCH,YAW RATE COMMANDS 

IF (TA LT. PIT(JX)) THEN 
JJX = JX - 1 

P1DDOT = (PIA(JX) - P1A(JJX))/(P1T(JX) - P1T(JJX)) 

P1DOT = PI DOT + PlDDOT*DT 
ELSE 

P1DDOT = (PIA(JX) - P1A(JJX))/(P1T(JX) - PIT(JJX)) 

P1DOT = P1DOT + PlDDOT*DT 

JX = JX+1 

ENDIF 

C WRITE(6,*) TA,P 1 DOT,P 1 DDOT 

IF (TA LT. P2T(JY)) THEN 
JJY = JY - 1 

P2DDOT = (P2A(JY) - P2A(JJY))/(P2T(JY) - P2T(JJY)) 

P2DOT = P2DOT + P2DDOT*DT 
ELSE 

P2DDOT = (P2A(JY) - P2A(JJY))/(P2T(JY) - P2T(JJY)) 

P2DOT = P2DOT + P2DDOT*DT 
JY = JY + 1 
END IF 

IF (TA LT. P3T(JZ)) THEN 
JJZ = JZ - 1 

P3DDOT = (P3A(JZ) - P3A(JJZ))/(P3T(JZ) - P3T(JJZ)) 

P3DOT = P3DOT + P3DDOT*DT 
ELSE 

P3DDOT = (P3A(JZ) - P3A(JJZ))/(P3T(JZ) - P3T(JJZ)) 

P3DOT = P3DOT + P3DDOT*DT 

JZ = JZ + 1 

ENDIF 

IF (TA .LT. P4T(JR)) THEN 
JJR = JR - 1 

P4DDOT = (P4A(JR) - P4A(JJR))/(P4T(JR) - P4T(JJR)) 
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P4DOT = P4DOT + P4DD0T*DT 
ELSE 

P4DD0T = (P4A(JR) - P4A(JJR))/(P4T(JR) - P4T(JJR)) 

P4DOT = P4DOT + P4DD0T*DT 

JR = JR + 1 

ENDIF 

IF (TA LT. P5T(JP» THEN 
J JP = JP - 1 

P5DD0T = (P5A(JP) - P5A(JJP))/(P5T(JP) - P5T(JJP)) 

P5D0T = P5D0T + P5DDOT*DT 
ELSE 

P5DDOT = (P5A(JP) - P5A(JJP))/(P5T(JP) - P5T(JJP)) 

P5D0T = P5D0T + P5DDOT*DT 

JP = JP + 1 

ENDIF 

IF (TA .LT. P6T(JYA)) THEN 
JJYA = JYA - 1 

P6DDOT = (P6A(JYA) - P6A(JJYA))/(P6T(JYA) - P6T(JJYA)) 

P6DOT = P6DOT + P6DDOT*DT 
ELSE 

P6DD0T = (P6A(JYA) - P6A(JJYA))/(P6T(JYA) - P6T(JJYA)) 

P6DOT = P6D0T + P6DD0T*DT 

JYA = JYA + 1 

ENDIF 

C ESTABLISH THE BODY POSITION, VELOCITY, AND ACCELERATION 
C OVER TIME IN INERTIAL REFERENCE FRAME COORDINATES 
C RELATIVE TO THE HOST 

PDOTO = VEX(1)*(1.- DEXP(-(VEX(2))*(TA**2))) 

PDOT(l) = PDOTO*DCOS(P(6)) + P1DOT 
PDOT(2) = PDOTO*DSIN(P(6)) + P2DOT 
PDOT(3) = PDOTO*DTAN(P(5)) + P3DOT 
PDOT(4) = P4DOT 
PDOT(5) = P5DOT 
PDOT(6) = P6DOT 

PDDOTO = VEX(1)*(2.*VEX(2)*TA*DEXP(-(VEX(2))*(TA**2))) 

PDDOT(l) = PDDOTO*DCOS(P(6)) + P1DDOT 

PDDOT(2) = PDDOTO*DSIN(P(6)) + P2DDOT 

PDDOT(3) = PDDOTO*DTAN(P(5)) + P3DDOT 

PDDOT(4) = P4DDOT 

PDDOT(5) = P5DDOT 
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PDDOT(6) = P6DD0T 

P(l) = P(l) + PDOT(l)*DT + PDD0T(1)*DT**2 
P(2) = P(2) + PDOT(2)*DT + PDDOT(2)*DT**2 
P(3) = P(3) + PD0T(3)*DT + PDDOT(3)*DT**2 
P(4) = P(4) + PD0T(4)*DT + PDDOT(4)*DT**2 
P(5) = P(5) + PD0T(5)*DT + PDDOT(5)*DT**2 
P(6) = P(6) + PD0T(6)*DT + PDDOT(6)*DT**2 

C CALCULATE COORDINATE TRANSFORM FROM INERTIA FRAME (HOST) 
C TO VEHICLE FRAME 

C TRANSPOSE OF THE SNAME 1950 CONVENTION (IAW ABKOWICZ) 

C ORDER IS YAW, PITCH, ROLL 

DS4 = DSIN(P(4)) 

DC4 = DCOS(P(4)) 

DS5 = DSIN(P(5)) 

DC5 = DCOS(P(5)) 

DS6 = DSIN(P(6)) 

DC6 = DCOS(P(6)) 

C 

CTRN(1,1) = DC5*DC6 
CTRN(1,2) = DC5*DS6 
CTRN(1,3) = -DS5 

CTRN(2,1) = -DS6*DC4 + DS4*DS5*DC6 
CTRN(2,2) = 004*006 + 054*085*086 
CTRN(2,3) = DS4*DC5 
CTRN(3,1) = DS4*DS6 +DC4*DC6*DS5 
CTRN(3,2) = -DS4*DC6 + DC4*DS5*DS6 
CTRN(3,3) = DC4*DC5 
C 

CTRN(4,4) = 1.0 
CTRN(4,5) = 0.0 
CTRN(4,6) = -DS5 
CTRN(5,4) = 0.0 
CTRN(5,5) = DC4 
CTRN(5,6) = DS4*DC5 
CTRN(6,4) = 0.0 
CTRN(6,5) = -DS4 
CTRN(6,6) = DC5*DC4 

C XFORM INERTIA RELATIVE VELOCITIES AND ACCELERATIONS 
C TO VEHICLE RELATIVE 
C CTRN (6X6) * PDOT (6X1) = E (6X1) 

C CTRN (6X6) *PDDOT (6X1) = EDOT (6X1) 
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C INITIALIZE THE E AND EDOT MATRICE TO ZERO AT EACH TIME STEP 



DO 110 1=1,6 
E(I) = 0.0 
EDOT(I) = 0.0 
110 CONTINUE 

DO 120 1=1,6 
DO 130 K=l,6 

E(I) = E(I) + PDOT(K)*CTRN(I,K) 

EDOT(I) = EDOT(I) + PDDOT(K)*CTRN(I,K) 

130 CONTINUE 
120 CONTINUE 

C CONVERT BODY VELOCITIES TO ADD IN HOST VESSEL SPEED EFFECTS 
C HOST SPEED IN (-) P(l) DIRECTION ONLY. ADD AS A POSITIVE SCALAR 

E(l) = E(l) + UHOST*DCOS(P(6))*DCOS(P(5)) 

E(2) = E(2) - UHOST*DSIN(P(6)) 

E(3) = E(3) + UHOST*DCOS(P(6))*DSIN(P(5)) 

C SOLVE FOR THE RIGHT HAND SIDE VECTOR WITH TOTAL FORCE 
C EDOT (6X1) * LHS (6X6) = RHS (6X1) 

C INITIALIZE THE RHSTOT TO ZERO AT EACH TIME STEP 
DO 140 1=1,6 
RHSTOT(I) = 0.0 
140 CONTINUE 

IF (JJ .EQ. 9) THEN 

WRITE(25,178) (acp(I),I=l,3) 

LABEL = 'CTRN' 

CALL MPRINT1(CTRN,6,6,9,LABEL) 

LABEL = 'LHS' 

CALL MPRINT 1(A,6,6,9,LABEL) 

ENDIF 
JJ = JJ + 1 

DO 200 1=1,6 
DO 220 K=l,6 

RHSTOT(I) = RHSTOT(I) +A(I,K)*EDOT(K) 
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220 CONTINUE 
200 CONTINUE 
AEDOT1 = RHSTOT(l) 

C CALCULATE THE RIGHT HAND SIDE OF THE HYDRODYNAMIC 
EQUATIONS 

C WITH NO CONTROL FORCE OR INTERACTION FORCE 
C SET CROSSFLOW SWITCH JJJ TO I 



JJJ = 1 

CALL HYDRO 

C CALCULATE THE INTERACTION FORCES IN FORCE MODEL ROUTINE 

CALL FM (FX,FY,FZ,MX,MY,MZ,AEDOTl,PTF,ATL,KODE) 

C CALCULATE THE CONTROL FORCES TO EFFECT THE 
C CHOSEN TRAJECTORY. ADD IN HYDRO FORCE 

RHSFM(l) = FX 
RHSFM(2) - FY 
RHSFM(3) = FZ 
RHSFM(4) = MX 
RHSFM(5) = MY 
RHSFM(6) = MZ 

DO 240 1=1,6 

RHSC(I) = RHSTOT(I) - RHSFM(I) 

240 CONTINUE 

C PRINT OUTPUT EVERY 10 DELTA T 



J = J+l 

IF (J .EQ. 10.) THEN 
J = 0 

WRITE(21,176) TA, (RHSTOT(I),I=l,6) 

WRITE(22,178) TA, KODE, (RHS(I), 1=1,6) 

WRITE(23,178) TA, KODE, (RHSFM(I),I=1,6) 

WRITE(24,177) TA, KODE, (RHSC(I),I=1,6) 

WRITE(25, 1 76) TA,(P(I),I= 1,3),(RTD*P(I),I=4,6) 

WRITE(26, 176) TA,(PDOT(I),I=l,3),(RTD*PDOT(I),I=4,6) 
WRITE(27,176) TA,(PDDOT(I),I=l,3),(RTD*PDDOT(I),I=4,6) 
WRITE(28, 176) TA,(E(I),I= 1 ,3 ),(RTD*E(I),I=4,6) 



70 



WRITE(29, 176) TA,(EDOT(I),I=l,3),(RTD*EDOT(I),I=4,6) 

176 F0RMAT(F5. 1, 1X,6(E1 1.4E3, IX)) 

177 F0RMAT(F6. 1, 1X,I2, 1X,6(E1 1.3E3, IX)) 

178 F0RMAT(F6. 1, IX, 12, 1X,6(E9.3, IX)) 

ENDIF 

RETURN 

END 
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£ ********************************************************************** 
SUBROUTINE FM (FX,FY,FZ,MX,MY,MZ, 

1 AEDOT 1 ,PTF, ATL,KODE) 



Q *********************** *** ******************** ***^/Qg***** ************ 

C SUBTOUTINE FORCE MODEL -- VEHICLE SPECIFIC FORCE AND MOMENT 
C MODELS ARE PLACED HERE 
C 

C CENTER OF X,Y,Z, COORDINATE SYSTEM IS ON HOST AXIS 
C IN THIS VERSION, IT IS PRESUMED THAT A NON-CRASH TRAJECTORY 
C SPECIFIED 

C MODEL IS VALID FOR THE STARBOARD HALF PLANE 
C (ALPHA CALCULATION) 

^ ******************************************** ****************** 



REAL* 8 MX,MY,MZ,MYC,MZC,MS,MTOT 

REAL* 8 P(6),E(6),XC, Y C,ZC,XMIN, YMIN,ZMIN 

REAL* 8 XMAX,YMAX,ZMAX,D2,HRAD,VRAD,UHOST,DIA 

REAL*8 RHOC,ATLMAX,VLEN,ACP(3),ACOF(6),BCP(3) 

REAL*8 BCOF(6),FX,FY,FZ,DOGS,PTF,PTFS,PTFC 
REAL* 8 UZ 1 LF,UZ 1 LM,UZ1 VF,UZ 1 VM,FACC,FACS,GAPS 
REAL*8 FYC,FZC,GAP,DOG,PT,SG,FTOT,POSRAD 
REAL* 8 UZ2LF,UZ2LM,X, Y,Z, ATL, ATLESTV, AXDRAG 
REAL* 8 U,V,W,PFQ,THETA,PSI,ALPFIA,FAREA,C1,C2,C4,C5 
REAL* 8 Dl, RAT, XS,YS,ZS, AEDOT 1,RHS1,TA 
REAL* 8 RHS(6) 

REAL* 8 DD( 1800) 

REAL*4 FOO,FOO 1 ,F002,PTF2 

INTEGER KODE 
COMMON /COMDD/ DD 

COMMON /FGEOM/ XC,YC,ZC,XMIN,XMAX,YMIN,YMAX,ZMIN,ZMAX, 
1 D2,HRAD, 

VRAD,UHOST,DIA,RHOC,ATLMAX,VLEN,ACP,ACOF,BCP,BCOF 
EQUIVALENCE (DD( 55), E(1)),(DD( 61), P(l)) 

EQUIVALENCE (DD( 10), TA),(DD( 71), RHS( 1)) 



PI = 3.14159 
DPR = 180.0/pi 
X = P(1) 

Y = P(2) 

Z = P(3) 

U = E(1) 

V = E(2) 
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W = E(3) 

PHI = P(4) 

THETA = P(5) 

PSI = P(6) 

ALPHA = ATAN2(Z,Y) 

FX = 0.0 
FY = 0.0 
FZ = 0.0 
MX = 0.0 
MY = 0.0 
MZ = 0.0 

FAREA = PI/4.*DIA**2 

Cl = ,5*RHOC*(UHOST**2)*FAREA 

C2 = CI*VLEN 

C4 = . 5 *(RHOC *U* *2)*FARE A 
C5 = C4*VLEN 
C6 = C5/12. 

C IN THIS SUBROUTINE, THE HOST IS TREATED AS A CYLINDER OF RADIUS 
HRAD. 

C POSRAD IS RADIUS FROM HOST AXIS TO SMALL VEHICLE ORIGIN. 

POSRAD = SQRT(Y*Y + Z*Z) 

D1 = SQRT( (X-XC)**2 + ( Y-YC)**2 + ( Z-ZC)**2 ) 

IF ((Z .LE. ZC) .AND. (Y .LE. YMAX) .AND. (Y .GE. YMIN) .AND. 

1 (X .LE. XMAX) .AND. (X .GE. XMIN)) THEN 



C VEHICLE IS IN THE LAUNCHWAY - ZONE 1 * 



ATL = ATLMAX-(SQRT( (X-XC)**2 + (Y-YC)**2 + (Z-ZC)**2))/VLEN 

ATLINV = ATL - ATLMAX 

AXDRAG= (.637341 +.853415*(ATLINV)**3)*C1 

IF (AXDRAG It. 0.0) AXDRAG - 0.0 

FTOT = AEDOT1 + AXDRAG 

PTF = ,5*(FTOT)/Cl 



FX = -AXDRAG 



UZ1LF = ((-2.04*ATL**3+5.59*ATL**2-4.49*ATL+1.13)*PTF + 
1 ( 1.40657-.55068*ATL)) 

FY = UZ1LF*C1 
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UZ1LM = -((. 153*ATL**3 198*ATL**2 - ,516*ATL +.508)*PTF) 

MZ = UZ1LM*C2 
KODE = 1 



ELSEEF ( (Z .LE. ZC) .AND. (Y .GE. YMAX) .AND. (X .GE. XMAX) 

1 .AND. (D1 .LT. D2)) THEN 

VEHICLE IS IN TRANSITION ZONE - HEMI- SPHERE ZONE 2 * 

* 



IF (D1 .LT. 0.001) D1 = 0.001 

RAT = D2/D1 

XS = XC + RAT*(X-XC) 

YS = YC + RAT*(Y-YC) 

ZS = ZC + RAT*(Z-ZC) 

FACC = (D2-D1)/D2 
FACS = D1/D2 

GAPS = SQRT(YS*YS + ZS*ZS) - HRAD - VRAD 

IF (GAPS .It. 0.001) GAPS = 0.001 

DOGS = DIA/GAPS 

ATL = ATLMAX 

ATLINV = 0.0 

AXDRAG = (.637341)*C1 

FXC = AEDOT1 + AXDRAG 

C AT LAUNCHER EXIT 
PTFC = ,5*(FXC)/C1 



FYC = UZ1LF*C1 

UZ1LF = ((-2.04*ATL**3+5.59*ATL**2-4.49*ATL+1. 13)*PTFC + 

1 (1.40657-.55068*ATL)) 

UZ1LM = -((. 153*ATL**3 -.198*ATL Ht *2 - .516^11. +.508)*PTFC) 
MZC = UZ1LM*C2 
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C PTF EFFECT AT SPHERE INTERSECTION POINT 



FXS = -RHS(l) + AEDOT1 

PTFS = ,5*(FXS)/C4 

UZ2LF = -(0.004*DOGS*PTFS) 

UZ2LM = (4. 13743 *DEXP(- 1 56. 1 97*( 1 ,/DOGS)))*PTFS 
Z2LMPTF = UZ2LM*C5 

C ADD WALL RELATIVE PITCH AND GAP EFFECT 

C APPLY EFFECTIVE PITCH ANGLE DUE TO LAUNCHWAY GEOMETRY 

PTEFF = ((-THETA* SIN( ALPHA) + PSI*COS( ALPHA)) * DPR) 

PT = PTEFF* (XS-XC)/(2.5*VLEN) 

CDBG WRITE(6,*) PTEFF, PT 

SG= SQRT(DOGS) 

UZ2LF = UZ2LF 

1 -((ACP( 1 )*DEXP(ACP(2)*( 1 ,/DOGS)**ACP(3)))+ 

1 (ACOF(l)*PT*SG+ACOF(2)*PT*SG*SG+ACOF(3)*PT*SG**3 + 

1 ACOF(4)*PT*PT*SG+ACOF(5)*PT*PT*SG*SG+ACOF(6)*PT*PT*SG**3)) 

FS = UZ2LF*C4 

UZ2LM = 

1 (BCOF(l)*PT*SG+BCOF(2)*PT*SG*SG+BCOF(3)*PT*SG**3 + 

1 BCOF(4)*PT*PT*SG+BCOF(5)*PT*PT*SG*SG+BCOF(6)*PT*PT*SG**3) 

MS = UZ2LM*C5+Z2LMPTF 

FX = FACC*(-FXC) + FACS*(-FXS) 

FY = FACC*FYC - FACS*FS*COS(ALPHA) + FACS*RHS(2) 

FZ = FACC*FZC - FAC S*FS*SIN( ALPHA) + FACS*RHS(3) 

MY = FACC*MYC - FAC S*MS*SIN( ALPHA) + FACS*RHS(5) 

MZ = FACC*MZC + FACS*MS*COS( ALPHA) + FACS*RHS(6) 

PTF = FACC*PTFC + FACS*PTFS 

KODE = 3 
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c 



ELSE 

£ ******************************************************************* 
C ZONE IS ZONE 2 * 



GAP = SQRT(Y*Y + Z*Z ) - HRAD - VRAD 
IF (GAP .LT. 0.001) GAP = 0.001 
DOG = DIA/GAP 
FX = -RHS(l) + AEDOT1 

C PTF EFFECT 

PTF - ,5*(-RHS(l) + AEDOTl)/C4 
UZ2LF = -(0.004*DOG*PTF) 

UZ2LM = (4.13743*DEXP(-156.197*(l./DOG)))*PTF 
Z2LMPTF = UZ2LM*C5 

C ADD WALL RELATIVE PITCH AND GAP EFFECT 

PT = (-THETA* SIN( ALPHA) + PSI*COS(ALPHA)) * DPR 
SG = SQRT(DOG) 

UZ2LF = UZ2LF 

1 -((ACP( 1)*DEXP(ACP(2)*( l./DOG)**ACP(3)))+ 

1 (ACOF(l)*PT*SG+ACOF(2)*PT*SG*SG+ACOF(3)*PT*SG**3 + 

1 ACOF(4)*PT*PT*SG+ACOF(5)*PT*PT*SG*SG+ACOF(6)*PT*PT*SG**3)) 

FTOT = UZ2LF*C4 

UZ2LM = 

1 (BCOF(l)*PT*SG+BCOF(2)*PT*SG*SG+BCOF(3)*PT*SG**3 + 

1 BCOF(4)*PT*PT*SG+BCOF(5)*PT*PT*SG*SG+BCOF(6)*PT*PT*SG**3) 

MTOT = UZ2LM*C5 + Z2LMPTF 

FX = -FX 

FY = - FTOT*COS(ALPHA) + RHS(2) 

FZ = - FTOT* SIN( ALPHA) + RHS(3) 
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MY = - MTOT*SIN( ALPHA) + RHS(5) 

MZ = MTOT*COS( ALPHA) + RHS(6) 

KODE = 2 

END IF 

RETURN 

END 
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£ He********************************************************************* 

SUBROUTINE HYDRO 

£ ********************************************************************** 
C DTRC 2510 HYDRODYNAMIC COEFFICIENTS 
C NOW INCLUDES ZERO FORWARD SPEED COEFFICIENTS 

C 

IMPLICIT NONE 

INTEGER* 2 JJJ,II,J 

INTEGER* 2 IDTAR 

REAL* 8 TEMP 1,TEMP2,TEMP3, TEMP 

REAL*8 UPFRAC,LOFRAC,URACER,URACERO,URACES,URACESO 
REAL*8 DRE3,DSPE3,DSSE3,DRE,DSPE,DSSE 
REAL* 8 DRIIS 1,DRHS2,DRHS3,DRHS4,DRHS5,DRHS6 
REAL* 8 ZRHS 1,ZRHS2,ZRHS3,ZRHS4,ZRHS5,ZRHS6 
REAL*8 BRHS1,BRHS2,BRHS3,BRHS4,BRHS5,BRHS6 
C 

REAL* 8 DD( 1800) 

COMMON /COMDD/ DD 
C 

REAL* 8 URUR,USUS,DELRCR,DELRCS,ACOF(8),BCOF(8),acp(4),bcp(4) 
REAL*8 DSP,DSS,DR,DB,UC,DT,RHS(6),E(6),P(6),QPROP,TA,YG,YB 
REAL* 8 PHISTERN, VATS, WATS, ETA,CETA,ETAM1 , BETAS 
REAL* 8 KBETA, WT,B OUY 

REAL* 8 ARAPPHI2,ARAPPHI4,CROSSLIM,SINTHETA 

REAL* 8 COSTHETA,SINPHI,COSPHI 

REAL* 8 UU,W,WW,PP,QQ,RR,PR,PQ,QR,UV,UW,UP 

REAL* 8 UQ,UR,VATTAIL,UUDB 

REAL *8 VP,VR,WP,WQ,SQRWWW,VSQRWWW 

REAL*8 WSQRWWW,CROSS(6),VORT(6) 

REAL* 8 DRDR,DSDS,UUDS,UUDR,UABW,UABV,WABW 
REAL*8 QABQ,VABV,RABR,FXP 
REAL* 8 PABP,DPDP,UUDP,UT, BETA, ALPHA, WR 
REAL* 8 ARRAY(6),WTINH20,DBDB 
C 

REAL* 8 XUDOT,XQQ,XRP,XVR,XWQ,XW,XDRDR,XDSSDSS 
REAL* 8 XDSPDSP,XDBDB 

REAL*8XW,XQ,XWDSP,XWDSS,XQDSS,XQDSP,XVDR,XRDR,XWW,XRR 

REAL* 8 XDRE2,XDSSE2,XDSPE2 

REAL* 8 YVDOT,YRDOT,YRU,YVU,YVABVU,YW 

REAL* 8 YDR,YDRN,YRABR,YVABV 

REAL* 8 YPU, YWP, YPABP, YPQ, Y STAR,YVAB VR,YHVAB V, YDRE, YDRE3 

REAL* 8 ZWDOT,ZQDOT,ZQU,ZWU,ZABWU,ZWW 

REAL* 8 ZDSS,ZDSSN,ZQABQ,ZWABW 

REAL* 8 ZVP,ZSTAR,ZPR,ZDSP,ZDSPN,ZHWABW,ZDB 

REAL* 8 ZDSSE,ZDSPE,ZDSPE3,ZDSSE3 
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REAL* 8 KPDOT,KVDOT,KRDOT,KQR,KPU,KRU,KWP,KVU,K4S 

REAL *8 K8S,KI,KDSP,KDSS 

REAL* 8 KPABP,KWR,KDR,KDRN,KSTAR 

REAL* 8 KDRE,KDRE3,KDSSE,KDSPE,KDSSE3,KDSPE3 

REAL* 8 MQDOT,MWDOT,MQU,MWU,MWABWR,MABWU,MWW,MDSS 

REAL* 8 MDSSN,MQABQ,MWABW 

REAL* 8 MRP,MSTAR,MVP,MDSP,MDSPN,MHWABW,MDB 
REAL* 8 MDSSE,MDSPE,MDSSE3,MDSPE3 

REAL* 8 NRDOT,NVDOT,NRU,NVU,NVABVR,NABVXJ,NW,NDR 
REAL* 8 NDRN,NVAB V,NRABR 

REAL* 8 NWP,NPQ,NSTAR,NPU,NHVABV,NDRE,NDRE3 
REAL* 8 IYY,IZZ,EXX,SLEN,MASS,ZB,ZG,XB,XG,XSTERN 
REAL*8XFORWARD,KTAINF,UAINF,IXZ,IXY,IYZ,YPDOT,NPDOT 
EQUIVALENCE (DD( 1), DT), (DD( 10), TA) 

EQUIVALENCE (DD( 11), JJJ), (DD( 16), BETA) 

EQUIVALENCE (DD( 15), ALPHA) 

EQUIVALENCE (DD( 17), FXP ),(DD( 1 8), QPROP) 

ETA ),(DD( 20), CETA) 

E(1)),(DD(61), P(l)) 

DSS),(DD( 68), DSP) 

DR),(DD( 71), RHS(l)) 

DB) 

EQUIVALENCE (DD(260),CROSS(1)),(DD(266), VORT(l)) 
EQUIVALENCE (DD(272),VATTAIL) 



EQUIVALENCE (DD( 19), 
EQUIVALENCE (DD( 55), 
EQUIVALENCE (DD( 67), 
EQUIVALENCE (DD( 69), 
EQUIVALENCE (DD(179), 



EQUIVALENCE (DD(301),XUDOT ),(DD(303),XW ) 

EQUIVALENCE (DD(305),XQ ),(DD(307),XQQ ) 

EQUIVALENCE (DD(309),XRR ),(DD(31 1),XRP ) 

EQUIVALENCE (DD(31 3 ),XVR ),(DD(315),XWQ ) 

EQUIVALENCE (DD(3 1 7),XW ),(DD(3 1 9),XWW ) 

EQUIVALENCE (DD(321),XDBDB ),(DD(323),XDRDR ) 
EQUIVALENCE (DD(325),XDSSDSS ),(DD(327),XDSPDSP ) 
EQUIVALENCE (DD(329),XDRE2 ),(DD(331),XDSSE2 ) 
EQUIVALENCE (DD(333),XDSPE2 ) 

EQUIVALENCE (DD(1703),XWDSP ),(DD(1705),XWDSS ) 
EQUIVALENCE (DD(1 707), XQDSP ),(DD(1709),XQDSS ) 
EQUIVALENCE (DD( 17 11), XVDR ),(DD(1713),XRDR ) 



EQUIVALENCE (DD(338),YVDOT ),(DD(340),YPDOT ) 
EQUIVALENCE (DD(342),YRDOT ),(DD(344),YPABP ) 
EQUIVALENCE (DD(346),YPQ ),(DD(348),YRU ) 

EQUIVALENCE (DD(350),YPU ),(DD(352),YVU ) 

EQUIVALENCE (DD(354),YVABVR ),(DD(356),YVABVU ) 
EQUIVALENCE (DD(358),YW ),(DD(360),YDR ) 

EQUIVALENCE (DD(362),YDRN ),(DD(364),YRABR ) 
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EQUIVALENCE (DD(366),YVABV ),(DD(368),YWP ) 

EQUIVALENCE (DD(3 70), YSTAR ),(DD(372),YHVABV ) 
EQUIVALENCE (DD(374),YDRE3 ),(DD(376),YDRE ) 

EQUIVALENCE (DD(3 90), ZWDOT ),(DD(392),ZQDOT ) 
EQUIVALENCE (DD(394),ZQU ),(DD(396),ZVP ) 

EQUIVALENCE (DD( 3 98), ZSTAR ),(DD(400),ZWU ) 

EQUIVALENCE (DD(402),ZABWU ),(DD(404),ZWW ) 

EQUIVALENCE (DD(406),ZDB ),(DD(408),ZDSS ) 

EQUIVALENCE (DD(410),ZDSSN ),(DD(412),ZDSP ) 

EQUIVALENCE (DD(414),ZDSPN ),(DD(416),ZQABQ ) 
EQUIVALENCE (DD(418),ZWABW ),(DD(420),ZHWABW ^ 
EQUIVALENCE (DD( 422),ZDSSE3 ),(DD( 424),ZDSPE3 ) 
EQUIVALENCE (DD( 426), ZDSSE ),(DD( 428),ZDSPE ) 
EQUIVALENCE (DD( 1 7 1 7),ZPR ) 

EQUIVALENCE (DD(440),KVDOT ),(DD(442),KPDOT ) 
EQUIVALENCE (DD(444),KRDOT ),(DD(446),KQR ) 

EQUIVALENCE (DD(448),KPABP ),(DD(450),KPU ) 

EQUIVALENCE (DD(452),KRU ),(DD(454),KWP ) 

EQUIVALENCE (DD(456),KSTAR ),(DD(458),KVU ) 

EQUIVALENCE (DD(460),KI ),(DD(462),K4S ) 

EQUIVALENCE (DD(464),K8S ),(DD(466),KDR ) 

EQUIVALENCE (DD(468),KDRN ),(DD(470),KDSS ) 

EQUIVALENCE (DD(472),KDSP ),(DD(474),KDRE ) 

EQUIVALENCE (DD(476),KDRE3 ),(DD(478),KDSSE ) 
EQUIVALENCE (DD(480),KDSPE ),(DD(482),KDSSE3 ) 
EQUIVALENCE (DD(484),KDSPE3 ) 

EQUIVALENCE (DD( 1 72 1 ),KWR ) 

EQUIVALENCE (DD(494),MWDOT ),(DD(496),MQDOT ) 
EQUIVALENCE (DD(498),MRP ),(DD(500),MQU ) 

EQUIVALENCE (DD(502),MWU ),(DD(504),MSTAR ) 

EQUIVALENCE (DD(506),MWABWR ),(DD(508),MABWU ) 
EQUIVALENCE (DD(510),MWW ),(DD(5 12),MDSS ) 

EQUIVALENCE (DD(514),MDSSN ),(DD(516),MDSP ) 
EQUIVALENCE (DD(518),MDSPN ),(DD(520),MDB ) 

EQUIVALENCE (DD(522),MQABQ ),(DD(524),MWABW ) 
EQUIVALENCE (DD( 526),MHWABW ),(DD( 528),MDSSE3 ) 
EQUIVALENCE (DD( 530),MDSPE3 ),(DD( 532),MDSSE ) 
EQUIVALENCE (DD( 534),MDSPE ) 

EQUIVALENCE (DD( 1 723 ),MVP ) 

EQUIVALENCE (DD(546),NVDOT ),(DD(548),NPDOT ) 
EQUIVALENCE (DD(550),NRDOT ),(DD(552),NPQ ) 



80 



non 



EQUIVALENCE (DD(554),NPU ),(DD(556),NRU ) 

EQUIVALENCE (DD(558),NSTAR ),(DD(560),NVU ) 

EQUIVALENCE (DD(562),NVABVR ),(DD(564),NABVU ) 
EQUIVALENCE (DD(566),NW ),(DD(568),NDR ) 

EQUIVALENCE (DD(570),NDRN ),(DD(572),NVABV ) 

EQUIVALENCE (DD(574),NRABR ),(DD(576),NWP ) 

EQUIVALENCE (DD(578),NFIVABV ),(DD(580),NDRE3 ) 
EQUIVALENCE (DD(582),NDRE ) 

EQUIVALENCE (DD(593),SLEN ),(DD(596),IYY ) 

EQUIVALENCE (DD(598),IZZ ),(DD(600),IXX ) 

EQUIVALENCE (DD(602),IXZ ),(DD(604),IYZ ) 

EQUIVALENCE (DD(606),IXY ),(DD(608),MASS ) 

EQUIVALENCE (DD(6 10), WT ),(DD(6 12),BOUY ) 
EQUIVALENCE (DD(6 14), ZB ),(DD(616),ZG ) 

EQUIVALENCE (DD(61 8), XB ),(DD(620),XG ) 

EQUIVALENCE (DD(622),XSTERN ),(DD(624),XFORWARD) 
EQUIVALENCE (DD(626),YG ),(DD(628),YB ) 

EQUIVALENCE (DD(629),WTINH20 ) 

EQUIVALENCE (DD(630),CROSSLIM) 

EQUIVALENCE (DD(631),LOFRAC ),(DD(632),UPFRAC ) 

EQUIVALENCE (DD(848),KTAINF ),(DD(849),UAINF ) 
EQUIVALENCE (DD(850),URACER0 ),(DD(851),URACES0 ) 
EQUIVALENCE (DD(858),URACER ),(DD(859),URACES ) 
EQUIVALENCE (DD(866),DRE ),(DD(867),DSSE ) 

EQUIVALENCE (DD(868),DSPE ) 

EQUIVALENCE (DD( 1 75 8), ARRAY( 1 )),(DD( 1 764), ACOF( 1 )) 
EQUIVALENCE (DD(1800),botdep), (DD(1772),BCOF(l)) 
EQUIVALENCE (DD(1780),acp(l)), (DD(1788), bcp(l)) 



FOR EFFICIENCY THE FOLLOWING VARIABLES ARE CALCULATED 

UU = E(1)*E(1) 

W = E(2)*E(2) 

WW = E(3)*E(3) 

UT = DSQRT(UU+VV+WW) 

PP = E(4)*E(4) 

QQ = E(5)*E(5) 

RR = E(6)*E(6) 

PR = E(4)*E(6) 

PQ = E(4)*E(5) 

QR = E(5)*E(6) 

UV = E(1)*E(2) 
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UW = E(1)*E(3) 

UP = E(1)*E(4) 

UQ = E(1)*E(5) 

UR = E(1)*E(6) 

VP = E(2)*E(4) 

VR = E(2)*E(6) 

WP = E(3)*E(4) 

WQ = E(3)*E(5) 

WR = E(3)*E(6) 

URUR = URACER* URACER 
USUS = URACES*URACES 
DELRCR = URACER - URACERO 
DELRCS = URACES - URACESO 
VATS = E(2) + XSTERN*E(6) 

WATS = E(3) - XSTERN*E(5) 

IF( DAB S(URACER) .GT. 0. 1) THEN 
DRE = DR - DATAN2( VATS, URACER ) 
ELSE 
DRE = DR 
ENDDF 

EF( DAB S(URACER) .GT. 0.1) THEN 
DSSE = DSS - DATAN2( WATS, URACES ) 
DSPE = DSP - DATAN2( WATS, URACES ) 
ELSE 

DSSE = DSS 
DSPE = DSP 
ENDDF 

DRE3 = DRE*DRE*DRE 

DSPE3 = DSPE*DSPE*DSPE 

DSSE3 = DSSE*DSSE*DSSE 

DRDR = DR* DR 

DBDB = DB*DB 

DSDS = DSS*DSS 

DPDP = DSP*DSP 

UUDR = UU*DR 

UUDB = UU*DB 

UUDP = UU*DSP 

UUDS = UU*DSS 

UABW = E(1)*ABS(E(3)) 

UABV = E(1)*ABS(E(2)) 

WABW = E(3)*ABS(E(3)) 

QABQ = E(5)*ABS(E(5)) 
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VABV = E(2)*ABS(E(2)) 

PABP = E(4)*ABS(E(4)) 

RABR = E(6)*ABS(E(6)) 

SQRWWW =DSQRT(W+WW) 

VSQRWWW = E(2) * SQRWWW 
WSQRWWW = E(3) * SQRWWW 
WTINH20 = WT - BOUY 
SINPIII = DSIN(P(4)) 

COSPHI = DC0S(P(4)) 

SINTHETA = DSIN(P(5)) 

COSTHETA = DC0S(P(5)) 

C 

IF( (DAB S( VATS) .GT. 0.001) .OR. (DABS(WATS) .GT. 0.001)) THEN 
TEMPI = VATS* VATS 
TEMP2 = WATS*WATS 
TEMP3 = TEMP 1+TEMP2 
ARAPPHI2 = TEMP3*2.*VATS*WATS/TEMP3 

ARAPPHI4 = TEMP3*4.*VATS*WATS*(TEMP2 - TEMP 1 )/(TEMP3 *TEMP3 ) 
BETAS = DATAN2(DSQRT(TEMP3),E(1)) 

PM1 STERN = DATAN2(WATS,-VATS) 

ELSE 

ARAPPHI2 = 0.0 
ARAPPHI4 = 0.0 
BETAS =0.0 
PH I STERN = 0.0 
END IF 
C 

IF((DABS(E(3)) .GT. 0.1D-10) .AND. (DABS(E(1)) .GT. 0.01)) THEN 
ALPHA = DATAN2(E(3),E( 1)) 

ELSE 

ALPHA = 0.0 
END IF 
C 

EF( (DABS(E(2)) .GT. 0.1D-20) .AND. (DABS(UT) .GT. .01)) THEN 
BETA = DASIN(-E(2)/UT) 

ELSE 

BETA = 0.0 
END IF 
C 

C FOR EFFICIENCY CALL THESE SUBROUTINES ONLY ONCE 
C 

IF ( JJJ EQ. 1) THEN 
CALL VORTEX 
CALL CROSFLOW 
CALL PROP 
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C 

DO 345 DDTAR = 1,10 
CALL T O WDARRY ( IDTAR) 

345 CONTINUE 
ENDIF 
C 

CETA 

ETAM1 = ETA - 1.0 

CALCULATE RIGHT HAND SIDES OF EQUATIONS OF MOTION 
AXIAL FORCES ZERO SPEED MODEL 
ZRHS1 = 0.0 

LATERAL FORCE ZERO SPEED MODEL 

ZRHS2 = YHVABV*VATS*DABS(VATS) 

NORMAL FORCE ZERO SPEED MODEL 

ZRHS3 - ZHWABW*WATS*DABS(WATS) 

ROLL MOMENT ZERO SPEED MODEL 

ZRHS4 = KPABP*PABP 

PITCH MOMENT ZERO SPEED MODEL 

ZRHS5 = MHWABW*WATS*DABS(WATS) 

YAW MOMENT ZERO SPEED MODEL 

ZRHS6 = NHVABV*VATS*DABS(VATS) 

CALCULATE DTRC HYDRODYNAMIC FORCES 
HYDRODYNAMIC AXIAL FORCES 

DRHS 1 = XQQ*QQ + XRR*RR + XRP*PR + XVR* VR + XWQ*WQ + XW* W 
1 + XW*UW + XWW*WW + XQ*UQ 

1 + (XDSSDSS*DSDS + XDSPDSP*DPDP + XDRDR*DRDR + 

1 + XDBDB*DBDB)*UU 

1 + (XWDSP*DSP + XWDSS*DSS)*UW + (XQDSP*DSP + XQDSS*DSS)*UQ 
1 + (XVDR*UV + XRDR*UR)*DR 
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HYDRODYNAMIC LATERAL FORCE 
DRHS2 = YPABP*PABP + YPQ*PQ 

1 + YRU*UR + YPU*UP + YVU*UV + YVABVU*UABV + YW*W + 

1 + YWP*WP 

1 + (YDR + YDRN*ETAM1)*UUDR + YVABV*VABV + YRABR*RABR + 

1 + YSTAR*UU 
1 + YVABVR*VSQRWWW 

HYDRODYNAMIC NORMAL FORCE 

DRHS3 = ZQU*UQ + ZWU*UW + ZABWU*UABW + ZVP*VP + ZSTAR*UU 
1 + (ZDSS + ZDSSN*ETAM1) * UUDS + (ZDSP + ZDSPN*ETAM 1 ) * UUDP 
1 + ZDB*UUDB 

1 + ZQABQ*QABQ + ZWABW*WABW + ZPR*PR + 

1 + ZWW*DAB S(E(3 ))* SQRWWW 

HYDRODYNAMIC ROLL MOMENT 

KBETA = BETAS*BETAS*(UU + VATS*VATS + WATS*WATS) 

DRHS4 = KPU*UP + KRU*UR + KQR*QR + KWP*WP 
1 + ((KDR + KDRN*ETAM1)*DR + KDSS*DSS + KDSP*DSP + KSTAR)*UU 
1 + KWR*WR + KVU*UV + KPABP*PABP + KI*VATTAIL + ARAPPffl4*K4S 
1 + KBETA*(K4S*DSIN(4.*PHISTERN) + K8S*DSIN(8.*PHISTERN)) 

HYDRODYNAMIC PITCH MOMENT 

DRHS5 = MQU*UQ + MQABQ*QABQ + MWABW*WABW 
1 + MVP* VP + MRP*PR + MSTAR*UU + MDB*UUDB 
1 + MWU*UW + MWAB WR* W SQRWWW + MABWU*UABW + 

1 +MWW*DABS(E(3))*SQRWWW 

1 + (MDSS + MDSSN*ETAM1)*UUDS + (MDSP + MDSPN*ETAM1)*UUDP 

HYDRODYNAMIC YAW MOMENT 

DRHS6 = NPQ*PQ + NPU*UP + NWP*WP + NSTAR*UU 
1 + NRU*UR + NVU*UV + NVABVR*VSQRWWW + NABVU*UAB V + 

1 + NW*W 

1 + (NDR + NDRN*ETAM1)*UUDR + NVABV*VABV + NRABR*RABR 

CALCULATE BODY FORCES, AND SUM FORCES INDEPENDENT OF ALPHA 
OR BETA 

C USED FOR BOTH LOW AND HIGH ANGLE OF ATTACK MODELS 
C 

C BODY AXIAL FORCE 
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BRHSl = -MASS*( WQ - VR - XG*(QQ + RR) + ZG*PR) 

1 - WTINH20*SINTHETA + FXP + ARRAY(l) 

1 + XDRE2*DRE*DRE*URUR 

1 + (XDSSE2*DSSE*DSSE + XDSPE2*DSPE*DSPE)*USUS 

BODY LATERAL FORCE 

BRHS2 = -MASS*(UR - WP + XG*PQ + ZG*QR ) 

1 + WTINH20*C0STHETA*SINPHI 
1 + CROSS(2) + V0RT(2) + ARRAY(2) 

1 + YDRE*( DELRCR*(-XSTERN*E(6) - E(2)) + DR*URUR) 

1 + YDRE3*DRE3*URUR 
1 + YDRE3*DRE3*URUR 

BODY NORMAL FORCE 

BRHS3 = -MASS*(VP - UQ + XG*PR - ZG*(PP + QQ)) 

1 + WTINH20*C0STHETA*C0SPHI 
1 + CROSS(3) + VORT(3) + ARRAY(3) 

1 + ZDSSE*( DELRCS*(XSTERN*E(5) - E(3)) + DSS*USUS) 

1 + ZDSPE*( DELRCS*(XSTERN*E(5) - E(3)) + DSP*USUS) 

1 + (ZDSPE3*DSPE3 + ZDSSE3*DSPE3)*USUS 

BODY ROLL MOMENT 

BRHS4 = MASS*ZG*(UR - WP) 

1 - (ZG*WT - ZB*BOUY)*COSTHETA*SINPHI 
1 + (YG*WT - YB*BOUY)*COSTHETA*COSPH3 
1 + QPROP + VORT(4) + ARRAY(4) 

1 + KDRE*( DELRCR*(XSTERN*E(6) - E(2)) + DR*URUR) 

1 + KDRE3*DRE3*URUR 

1 + KDSSE * (DELRCS*(XSTERN*E(5) - E(3)) + DSS*USUS) 

1 + KDSPE * (DELRC S *(X STERN *E(5 ) - E(3)) + DSP*USUS) 

1 + (KDSPE3*DSPE3 + KDSSE3*DSPE3)*USUS 

BODY PITCH MOMENT 

BRHS5 = (-IXX + IZZ)*PR + MASS*(-ZG*(WQ - VR) + XG*(VP - UQ)) 
1 - (XG*WT - XB*BOUY)*COSTHETA*COSPHI 
1 - (ZG*WT - ZB*BOUY)*SINTHETA 
1 + CROSS(5) + VORT(5) + ARRAY(5) 

1 + MDSSE * (DELRCS*(XSTERN*E(5) - E(3)) + DSS*USUS) 

1 + MDSPE * (DELRCS*(XSTERN*E(5) - E(3)) + DSP*USUS) 

1 + (MDSPE3*DSPE3 + MDSSE3*DSPE3)*USUS 
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BODY YAW MOMENT 

BRHS6 = (-IYY + DCX)*PQ - MASS*(XG*(UR - WP)) * 

1 + (XG*WT - XB*BOUY)*COSTHETA*SINPHI 
1 + (YG*WT - YB*BOUY)*SINTFEETA 
1 + (XG*WT - XB*(BOUY))*COSTHETA*SINPHI 
1 + (YG*WT - YB*(BOUY))*SINTHETA 
1 + CROSS(6) + V0RT(6) + ARRAY(6) 

1 + NDRE*(DELRCR*(-XSTERN*E(6) - E(2)) + DR*URUR) 

1 + NDRE3*DRE3*URUR 

SUM ALL FORCES 

RHS(l) = BRHS1 + LOFRAC*DRHS 1 + UPFRAC*ZRHS1 
RHS(2) = BRHS2 + LOFRAC*DRHS2 + UPFRAC*ZRHS2 
RHS(3) = BRHS3 + LOFRAC*DRHS3 + UPFRAC*ZRHS3 
RHS(4) = BRHS4 + LOFRAC*DRHS4 + UPFRAC*ZRHS4 
RHS(5) = BRHS5 + LOFRAC *DRHS 5 + UPFRAC*ZRHS5 
RHS(6) = BRHS6 + LOFRAC *DRHS6 + UPFRAC*ZRIIS6 



RETURN 

END 
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Appendix B. Trajectory input files. 



This appendix presents the geometry limitation file and trajectory command files used to 
generate the simulation trajectories. 
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Geometry Limitation File GEOM.LIM 



93.57 16.50 0. XC YC ZC [FT] 

85.275 15.04 0. XMIN YMIN ZMIN [FT] 

93.576 16.50 0. XMAX YMAX ZMAX [FT] 

8.00 16.5 0.8125 10.5 D2 HRAD VRAD SLEN [FT] 

0.1 1.625 2.0 1.47 VHOST [FT/SEC] DIA [FT] RHOC [SLUG/FT A 3] 

ATLMAX 
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Traj 1 



2.0 0.04 VELOCITY AND ACCELERATION EXPONENTS 

XVEL,YVEL,ZVEL CORRECTION COMMANDS [FT/SEC] [SEC] 



P1A 


PIT 


P2A 


P2T 


P3A 


P3T 


0.00 


0.00 


0.00 


0.00 


0.00 


0.00 


0.00 


20.00 


0.00 


600.00 


0.00 


600.00 


0.00 


30.00 


0.00 


700.00 


0.00 


700.00 


0.00 


35.00 


0.00 


750.00 


0.00 


720.00 


0.00 


40.00 


0.00 


800.00 


0.00 


750.00 


0.00 


799.00 


0.00 


850.00 


0.00 


800.00 


0.00 


899.00 


0.00 


860.00 


0.00 


810.00 


0.00 


900.00 


0.00 


870.00 


0.00 


820.00 


0.00 


990.00 


0.00 


880.00 


0.00 


850.00 


0.00 


999.00 


0.00 


900.00 


0.00 


900.00 



ROLL, PITCH, YAW RATE COMMANDS [DEG/SEC] [SEC] 



P4A 


P4T 


P5A 


P5T 


P6A 


P6T 


0.00 


0.00 


0.00 


0.00 1 


0.00 0.00 


0.00 


600.00 


0.00 


600.00 


0.00 


20.00 


0.00 


700.00 


0.00 


700.00 


0.00 


25.00 


0.00 


750.00 


0.00 


750.00 


0.00 


28.00 


0.00 


800.00 


0.00 


800.00 


0.00 


750.00 


0.00 


900.00 


0.00 


850.00 


0.00 


800.00 


0.00 


999.00 


0.00 


860.00 


0.00 


810.00 


0.00 


999.00 


0.00 


870.00 


0.00 


820.00 


0.00 


999.00 


0.00 


880.00 


0.00 


850.00 


0.00 


999.00 


0.00 


900.00 


0.00 


900.00 
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Traj_2 



2.0 0.04 VELOCITY AND ACCELERATION EXPONENTS 

XVEL,YVEL,ZVEL CORRECTION COMMANDS [FT/SEC] [SEC] 



P1A 


PIT 


P2A 


P2T 


P3A 


P3T 


0.00 


0.00 


0.00 


0.00 


0.00 


0.00 


0.00 


20.00 


0.00 


600.00 


0.00 


600.00 


0.00 


30.00 


0.00 


700.00 


0.00 


700.00 


0.00 


35.00 


0.00 


750.00 


0.00 


720.00 


0.00 


40.00 


0.00 


800.00 


0.00 


750.00 


0.00 


799.00 


0.00 


850.00 


0.00 


800.00 


0.00 


899.00 


0.00 


860.00 


0.00 


810.00 


0.00 


900.00 


0.00 


870.00 


0.00 


820.00 


0.00 


990.00 


0.00 


880.00 


0.00 


850.00 


0.00 


999.00 


0.00 


900.00 


0.00 


900.00 



ROLL, PITCH, YAW RATE COMMANDS [DEG/SEC] [SEC] 



P4A 


P4T 


P5A 


P5T 


P6A 


P6T 


0.00 


0.00 


0.00 


0.00 1 


0.00 0.00 


0.00 


600.00 


0.00 


600.00 


0.00 


20.00 


0.00 


700.00 


0.00 


700.00 


-2.00 


25.00 


0.00 


750.00 


0.00 


750.00 


0.00 


30.00 


0.00 


800.00 


0.00 


800.00 


0.00 


35.00 


0.00 


900.00 


0.00 


850.00 


0.00 


36.00 


0.00 


999.00 


0.00 


860.00 


0.00 


810.00 


0.00 


999.00 


0.00 


870.00 


0.00 


820.00 


0.00 


999.00 


0.00 


880.00 


0.00 


850.00 


0.00 


999.00 


0.00 


900.00 


0.00 


900.00 
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Traj_3 



2.0 0.04 VELOCITY AND ACCELERATION EXPONENTS 
XVEL,YVEL,ZVEL CORRECTION COMMANDS [FT/SEC] [SEC] 



P1A 


PIT 


P2A 


P2T 


P3A 


P3T 


0.00 


0.00 


0.00 


0.00 0.00 1 


0.00 


0.00 


10.00 


0.00 


600.00 


0.00 


600.00 


0.00 


15.00 


0.00 


700.00 


0.00 


700.00 


0.00 


20.00 


0.00 


750.00 


0.00 


720.00 


-4.00 


30.00 


0.00 


800.00 


0.00 


750.00 


-4.00 


799.00 


0.00 


850.00 


0.00 


800.00 


-4.00 


899.00 


0.00 


860.00 


0.00 


810.00 


-4.00 


900.00 


0.00 


870.00 


0.00 


820.00 


-4.00 


990.00 


0.00 


880.00 


0.00 


850.00 


-4.00 


999.00 


0.00 


900.00 


0.00 


900.00 



ROLL,PITCH,YAW RATE COMMANDS [DEG/SEC] [SEC] 



P4A 


P4T 


P5A 


P5T 


P6A 


P6T 


0.00 


0.00 


0.00 


0.00 i 


0.00 0.00 


0.00 


600.00 


0.00 


600.00 


0.00 


20.00 


0.00 


700.00 


0.00 


700.00 


0.00 


25.00 


0.00 


750.00 


0.00 


750.00 


0.00 


28.00 


0.00 


800.00 


0.00 


800.00 


0.00 


750.00 


0.00 


900.00 


0.00 


850.00 


0.00 


800.00 


0.00 


999.00 


0.00 


860.00 


0.00 


810.00 


0.00 


999.00 


0.00 


870.00 


0.00 


820.00 


0.00 


999.00 


0.00 


880.00 


0.00 


850.00 


0.00 


999.00 


0.00 


900.00 


0.00 


900.00 
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Traj_4 

Traj_4 uses Trajl with Vhost = 0.01 in Geom.lim 
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Appendix C. Model Fit Analysis of Variance. 
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A/L 


Variance 

Source 


Degrees of 
Freedom 


Sum of 
Squares 


Mean Square 


F Ratio 


Prob > F 


.68 


Model 


1 


0.04653676 


0.046537 


4.0022 


0.0627 




Error 


16 


0.18604585 


0.011628 








Total 


17 


0.23258261 








.84 


Model 


1 


0.13460215 


0.134602 


2.6603 


0.1339 j 




Error 


16 


0.50595732 


0.050596 








Total 


11 


0.64055947 








1.0 


Model 


1 


11.940734 


11.9407 


25.7911 


<0.0001 




Error 


70 


32.408566 


0.4630 








Total 


71 


44.349300 








1.16 


Model 


1 


3.2888807 


3.28888 


88.9794 


<0.0001 




Error 


16 


0.5913961 


0.03696 








Total 


17 


3.8802768 








1.32 


Model 


1 


4.7664514 


4.76645 


166.9032 


<0.0001 




Error 


34 


0.9709780 


0.02856 








Total 


35 


5.9374294 








1.46 


Model 


1 


0.46515718 


0.465157 


31.1761 


<0.0002 




Error 


10 


0.14920309 


0.01492 








Total 


11 


0.61436 









Table C-l ANOVA Table for PTF Effect on Launchway Force 



A/L 


Variance 

Source 


Degrees of 
Freedom 


Sum of 
Squares 


Mean Square 


F Ratio 


Prob > F 


.68 


Model 


1 


0.04653676 


0.046537 


4.0022 


0.0627 




Error 


16 


0.18604585 


0.011628 








Total 


17 


0.23258261 








.84 


Model 


1 


0.13460215 


0.134602 


2.6603 


0.1339 




Error 


16 


0.50595732 


0.050596 








Total 


11 


0.64055947 








1.0 


Model 


1 


11.940734 


11.9407 


25.7911 


<0.0001 




Error 


70 


32.408566 


0.4630 








Total 


71 


44.349300 








1.16 


Model 


1 


3.2888807 


3.28888 


88.9794 


<0.0001 




Error 


16 


0.5913961 


0.03696 








Total 


17 


3.8802768 








1.32 


Model 


1 


4.7664514 


4.76645 


166.9032 


<0.0001 




Error 


34 


0.9709780 


0.02856 








Total 


35 


5.9374294 








1.46 


Model 


1 


0.46515718 


0.465157 


31.1761 


<0.0002 




Error 


10 


0.14920309 


0.01492 








Total 


11 


0.61436 









Table C-2 ANOVA Table for PTF effect on Launchway Moment 
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